Example of CalibrationΒΆ
Import stuff
import warnings
import pandas as pd
import matplotlib.pyplot as plt
from dask.distributed import Client
from dask_jobqueue import SLURMCluster as Cluster
import nevergrad as ng
from mhpc_project.matsch_b2 import CalibrationModel, Variables, Loss, NGO
from geotopy.utils import date_parser, comparison_plot, ProgressBar, ParametersLogger
from geotopy.measures import KGE
Launch Dask cluster using dask-jobqueue
cluster = Cluster()
cluster.adapt(maximum_jobs=128)
client = Client(cluster)
print(cluster.job_script())
#!/usr/bin/env bash
#SBATCH -J dask-worker
#SBATCH -e /dev/shm/dask-worker-%J.err
#SBATCH -o /dev/shm/dask-worker-%J.out
#SBATCH -p regular2
#SBATCH -n 1
#SBATCH --cpus-per-task=4
#SBATCH --mem=2GB
#SBATCH -t 00:20:00
#SBATCH --hint=nomultithread
#SBATCH --hint=compute_bound
export NUMEXPR_NUM_THREADS=1
export NUMEXPR_MAX_THREADS=1
export OMP_NUM_THREADS=1
export TMPDIR=/dev/shm
/home/scampane/scratch/MHPC-project/.env/bin/python -m distributed.cli.dask_worker tcp://10.6.0.4:40470 --nthreads 4 --memory-limit 2.00GB --name name --nanny --death-timeout 60 --local-directory $TMPDIR --lifetime 18m --lifetime-stagger 1m --interface ib0
Load observations, model, loss function and calibration strategy
observations = pd.read_csv('../data/Matsch B2/obs.csv',
na_values=['-9999', '-99.99'],
usecols=[0, 7],
parse_dates=[0],
date_parser=date_parser,
index_col=0,
squeeze=True)
observations.index.rename('datetime', inplace=True)
model = CalibrationModel('../data/Matsch B2/geotop', run_args={'timeout': 120})
variables = Variables('../data/Matsch B2/variables.csv')
measure = KGE(observations)
loss = Loss(model, variables, measure)
calibration = NGO(loss, budget=4096, num_workers=512)
Evaluate the model and plot the risults
simulation = model()
print(f"Before optimization loss is {measure(simulation)}")
comparison_plot(observations,
simulation,
desc='Soil moisture content @ 5cm',
rel=True)
plt.show()
Before optimization loss is 0.7523503598701219
Run calibration
logger = ParametersLogger(loss.massage)
calibration.optimizer.register_callback('tell', logger)
calibration.optimizer.register_callback('tell', ProgressBar(total=calibration.optimizer.budget))
with warnings.catch_warnings():
warnings.simplefilter("ignore")
recommendation, best_loss = calibration(executor=client)
(256_w,512)-aCMA-ES (mu_w=131.7,w_1=2%) in dimension 21 (seed=<module 'time' (built-in)>, Mon Nov 23 22:06:31 2020)
Plot results and show table of calibrated values
simulation = model(**recommendation)
print(f"After optimization loss is {best_loss}")
comparison_plot(observations,
simulation,
desc='Soil moisture content @ 5cm',
rel=True)
plt.show()
After optimization loss is 0.49338997022896586
pd.DataFrame.from_dict(recommendation, orient='index', columns=['value'])
| value | |
|---|---|
| NVanGenuchten | 1.293041 |
| AlphaVanGenuchten | 0.000552 |
| ThetaSat | 0.578172 |
| ThetaRes | 0.066367 |
| NormalHydrConductivity | 0.000276 |
| SoilRoughness | 168.386071 |
| LSAI | 4.317711 |
| CanopyFraction | 0.289856 |
| DecayCoeffCanopy | 16.953715 |
| RootDepth | 286.314931 |
| MinStomatalRes | 67.178791 |
| VegReflectVis | 0.136903 |
| VegReflNIR | 0.295314 |
| VegTransVis | 0.017572 |
| VegTransNIR | 0.089866 |
| CanDensSurface | 0.932729 |
| SoilAlbVisDry | 0.063979 |
| SoilAlbNIRDry | 0.329262 |
| SoilAlbVisWet | 0.215956 |
| SoilAlbNIRWet | 0.305602 |
| SoilEmissiv | 0.929224 |
logger.parallel_coordinate_plot()
Loading HiPlot...